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Abstract 

We present results of an ongoing study of the nature of the high temperature 
crossover in QCD with two light fermion flavors. These results are obtained 
with the conventional staggered fermion action at the smallest lattice spacing 
to date — approximately 0.1 fm. Of particular interest are (1) a study of the 
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temperature of the crossover, an important indicator of continuum scaling, 
(2) a determination of the induced baryon charge and baryon susceptibility, 
used to study the dissolution of hadrons at the crossover, (3) the scalar sus- 
ceptibility, a signal for the appearance of soft modes, and (4) the chiral order 
parameter, used to test models of critical behavior associated with chiral sym- 
metry restoration. From our new data and published results for N t = 4, 6, 
and 8, we determine the QCD magnetic equation of state from the chiral order 
parameter using 0(4) and mean field critical exponents and compare it with 
the corresponding equation of state obtained from an 0(4) spin model and 
mean field theory. We also present a scaling analysis of the Polyakov loop, 

suggesting a temperature dependent "constituent quark free energy." 
12.38.Gc,11.15.Ha,12.38.Mh,12.38.Aw,24.85.+p,64.30+t 
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I. INTRODUCTION 



Lattice simulations of high temperature QCD provide, at present, our only firmly 
grounded theoretical insights into the phenomenology of the transition from hadronic mat- 
ter to the quark-gluon plasma, and into the nature of the plasma itself. Our ultimate goals 
in lattice simulations at high temperature include (1) establishing continuum scaling, im- 
portant not only for determining the temperature of the crossover, but necessary for the 
validity of all dynamical fermion simulations, (2) establishing the mechanism for the dis- 
solution/formation of hadrons at the crossover, (3) exploring the intricate critical behavior 
associated with the phase transition, and (4) obtaining a quantitative characterization of 
the quark-gluon plasma, including the equation of state. To achieve these goals requires a 
combination of advances in algorithms and computing power . Recent improvements in 
quenched lattice algorithms hold promise for dynamical fermion simulations HQ. Here we 
present results of simulations with the conventional staggered fermion action at the smallest 
lattice spacing to date. 

The most detailed previous simulations with two-flavors of staggered fermions were car- 
ried out at N t = 4 (a « 0.3 fm) [|, N t = 6 (a « 0.2 fm) and |JV t = 8(a« 0.15 fm) 
0. Results of those simulations reaffirmed the hypothesis that for two flavors there is a 
rapid crossover, but no phase transition at nonzero quark mass. Moreover, the ratio T c /m p 
was found to be consistent with previous measurements. Despite the lack of surprises at 
N t = 6 and 8, there are strong reasons to push to still smaller lattice spacing. The full 
flavor symmetry in the staggered fermion scheme is restored only in the continuum limit. 
Sigma models of chiral symmetry restoration suggest that the crossover becomes a first order 
phase transition as the number of flavors is increased || . If flavor symmetry breaking were 
to cause an effective undercounting of quark flavors, one might expect a more pronounced 
crossover, or even a genuine phase transition at smaller lattice spacing. Moreover, strong 
coupling distortions in the hadron spectrum at Nt = 6 undermine credibility in a determi- 
nation of the crossover temperature or even in the plausibility of having achieved a scaling 
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ratio T c jm p . For these reasons we undertook a simulation at a ~ 0.1 fm. Preliminary results 
were reported at the Bielefeld and Melbourne Lattice conferences P,PJT0|j. 



The simulation with two flavors was carried out using the R algorithm described in (IT 
at two quark masses, am q = 0.008 and 0.016 and six couplings, 6/g 2 = 5.65, 5.70, 5.725, 
5.75, 5.80, 5.85, except that the 5.725 coupling was not simulated at the higher mass. In each 
case the simulation was extended to at least 2000 molecular dynamics time units. Lattices 
were saved at intervals of 8 time units. For the present analysis the first 500 time units were 
omitted for equilibration. Most of the results reported here are based on an analysis of the 
approximately 180 remaining lattices at each parameter pair. For the sake of comparison, 
we also present some new results for N t = 6, based on a compilation of lattices from our 
equation of state study Jl2| . 



In addition to the standard variables, the chiral order parameter and the Polyakov loop, 
we introduce a fuzzy Polyakov loop in an effort to regulate the ultraviolet divergence that 
becomes increasingly troublesome at small lattice spacing. We also measure the baryon 
susceptibility, the induced quark number, and the disconnected chiral susceptibility. Results 
are presented in Sec. 2. We show that with our choice of masses and couplings, the crossover 
is more clearly defined in terms of the baryon susceptibility and induced quark number. In 
Sec. 3 we present a scaling analysis of the Polyakov loop variable for world data ranging 
from N t — 4 to N t — 12, in terms of a constituent quark free energy. We also present a 
critical scaling analysis of the dependence of (ipipj on temperature and quark mass (magnetic 
equation of state), including a compilation of world data for this variable. This analysis is 
carried out in the context of mean field as well as 0(4) critical behavior. 



II. LOCATING THE CROSSOVER 
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A. Polyakov loop, (ipip), fuzzy loop 



The high temperature crossover is conventionally located from the inflection point in a 
plot of the Polyakov loop or chiral order parameter as a function of coupling. Figures [l] 
and |2| plot these quantities at fixed bare quark mass as a function of coupling. It is clear 
that despite the small errors, at these quark masses the effect of the the crossover on these 
quantities is subtle, indeed. 

The Polyakov loop measures the change f(T,m q ) in the ensemble free energy caused by 
the introduction of a static test quark. In particular, with our normalization 

(ReP) = 3exp[-/(T,m (? )/T] (2.1) 

As the lattice spacing shrinks, the test quark self-energy develops an ultraviolet divergence, 
which may overwhelm a crossover signal in this quantity. Thus one might hope that by 
increasing the radius of the test charge, one might regulate the divergence and recover a 
signal at least as strong as was seen on a coarser lattice. Accordingly, we constructed a 
"fuzzy Polyakov loop" variable in analogy with techniques introduced for glueball sources 
|T3| . To do so, we replaced the conventional product of forward links 

(HeP(x)> = /*& nV t (x,*)\ (2-2) 

with 

(ReF(x)> = /Tr U ftfat)) ( 2 -3) 

where 

F t (x,t) = aU t + f3j2U st ^ e (2.4) 

The six staples associated with the link £7t(x, t) are the usual three- link products of the form 

C4(x, t)U t {x + x, t)f/J(x, t + 1) (2.5) 

centered on the link from (x, t) to (x, t + 1). A weighting a = (3 = 1/7 was determined from 
a rough optimization of the variance. By construction the weights sum to 1 and we do not 



project the matrix F t onto SU(3), so this observable still creates a source of precisely one 
color triplet. 



that smearing indeed reduces the relative statistical error, but apparently does not enhance 
the crossover signal. 



The conventional Feynman path integral simulates the grand canonical ensemble in 
baryon number at zero chemical potential. The baryon susceptibility measures fluctuations 
in the baryon number of the ensemble. It is defined as the derivative of the baryon charge 
density with respect to chemical potential. The susceptibility can be defined separately for 
each flavor. Thus with two quark flavors, two susceptibilities can be measured: a flavor 
singlet and flavor nonsinglet fl5(| . Both quantities are obtained at zero chemical potential. 



The nonsinglet susceptibility is compared with results for lower N t [[HI in Fig. f|. Also 
indicated is the free lattice quark value for each N t . A common feature is the abrupt rise 
in susceptibility at crossover, followed by an asymptotic approach to the free lattice quark 
value. At lower values of Nt, where the crossover has been located with traditional methods, 
we find that the baryon susceptibility reaches 1/3-1/2 of the free quark asymptotic value at 
the crossover. Since this observable is based on a conserved charge, it is not renormalized, 
and offers our most distinctive signal for the crossover. Following the same rule for N t = 12 
places the crossover in the am q = 0.008 series between Q/g 2 = 5.65 and 5.70 and in the 
am q = 0.016 series between 6/g 2 = 5.75 and 5.80. 

To convert these results to a temperature, we use a scale in which the rho meson mass 
is taken to be 770 MeV, regardless of quark mass. The rho mass in lattice units is obtained 
in turn from an empirical fit to masses in a compilation of two-flavor staggered fermion 



Figure ^ gives the result. Comparing with the standard Polyakov loop in Fig. ^|, we see 



B. Baryon susceptibility 
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spectral simulations |f24}| , and includes an extrapolation beyond the parameter range of 
spectral measurements (Q/g 2 > 5.7) using tadpole-improved asymptotic scaling. Details 
are given in the Appendix. From the baryon susceptibility we then place the crossover at 
T c = 143 — 154 MeVp at the lighter quark mass and T c = 142 — 150 MeV p at the heavier 
quark mass. 



C. Induced quark number 



Another confinement-sensitive observable is the induced quark number [|17]]. This ob- 
servable measures the total residual light-quark number in an ensemble containing a single 
test quark. 

Qind = J Pind(r)d 3 r (2.7) 

The induced quark number density p in d is measured in the presence of the test quark. 
Operationally, the quantity measures the correlation between the Polyakov loop and the 
light-quark density |p"7f . 

This observable is subject to considerable fluctuation. We found it particularly effective 
to introduce the test charge through the fuzzy Polyakov loop variable described above. The 
baryon density of the dynamical quark is computed using a random source estimator |17| . 



To improve the signal further, we adjusted the number of random sources using an adaptive 
procedure on a configuration-by-configuration basis as follows: starting with a minimum of 
20 random sources the variance of the total induced baryon number was estimated. If the 
variance was greater than tolerance, another 20 random sources were added to the sample, 
and so on, up to a maximum of 80 for one configuration. 

Results are shown in Fig. [5] and compared with results from simulations at lower N t . The 
induced quark number is expected to be exactly —1 at zero temperature, since confinement 
requires screening of the test charge by a single antiquark. At the crossover, this quantity 
rises rapidly, approaching zero in the high temperature phase. At lower N t the induced quark 
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number reaches approximately —0.1 at crossover. Thus the induced quark number density 
can give an operational definition of the crossover. Applying this rule to the N t = 12 
data gives Q/g 2 = 5.65 — 5.70 (T = 143 — 154 MeV p ) in the am q = 0.008 series and 
Q/g 2 = 5.70 — 5.80 (T = 134 — 150 MeV p ) in the am q = 0.016 series. This crossover location 
is consistent with, but somewhat less precise than, that found from the baryon susceptibility. 



D. Chiral Susceptibility 

Another signal for the crossover is the singlet chiral susceptibility 



Xr, 



dm 



(2.8) 

6/9 2 



which measures fluctuations in the chiral order parameter (ipip) ||. Here we differentiate 
with respect to equated up and down quark masses. Since the derivative is the space-time 
integral of the correlator ^^(^^^(O)^, a peak in this observable occurs at a minimum 
in the a meson screening mass, indicating the presence of a soft mode. Such soft modes 
are expected in models of critical behavior ||. Like the a meson propagator, the singlet 
chiral susceptibility can be decomposed into two contributions: quark-line-connected and 
quark-line-disconnected: 

Xm Xconn Xdiso (2.9) 



The connected susceptibility in another guise is the derivative of (ipip) with respect to va- 
lence quark mass at fixed sea quark mass and Q/g 2 ||18|| , while the disconnected susceptibility 
is the configuration variance of (tpip)'- 



Xdi 



ISC 



u 

where {jpipj f is the expectation value computed on a single gauge configuration and () 
denotes averaging over all gauge configurations in the sample. 

We used the standard Gaussian random source estimator for (ipip) 



One must take care in using this estimator for Eq. ( [2.10| ), since the square involves cor- 
related quartic terms in the random source variable. Averaging then reintroduces part of 
the connected susceptibility along with the disconnected susceptibility. We chose to avoid 
the problem altogether by computing the estimator from five separate random sources on 
each configuration and forming products in Eq. (|2.10|) only from pairs of estimators with 
different random sources. (For the smaller volume N t = 6 data we used 33 random sources 
per configuration.) 

Results for x disc and Xconn are compared in Figs. || and [7| with the results of Karsch and 
Laermann at N t = 4 || and results from a reanalysis of a set of Nt = 6 lattices retained 
from an equation of state study |Ll|. The benefits of the larger sample in the N t = 4 
data are clearly evident in this variable. At the lighter quark mass (am q /T = 0.08) the 
Nt — 4 octagons define a crossover peak at about 145 MeV p . This peak can be sought in 
the N t = 6 data at a comparable quark mass (am q /T = 0.075). The data are crude but 
not inconsistent with a comparable peak. The bursts plot N t = 12 data at the lighter quark 
mass (am q /T = 0.096), which are not inconsistent with a peak in the range 140-160 MeV p . 
The signal in the higher mass Nt = 12 data (diamonds) is inconclusive. 

That the signal in Xdisc should be weaker at higher quark mass is expected from models 
of critical behavior that place the critical point at zero quark mass. The mass ratio m^/m 2 p 
measures proximity to the critical point. In the N t = 4 data, this ratio is less than 0.1 for 
the lighter quark mass (am q = 0.020) data and between 0.1 and 0.2 for the am q = 0.0375 
data. For the N t = 6 data it is also between 0.1 and 0.2. For the N t = 12 data, it is between 
0.2 and 0.3 for am q = 0.008 and 0.3 to 0.4 for am q = 0.016. Thus we expect a weakening of 
the crossover signal at the quark masses treated in our higher Nt data. 

The connected part of the susceptibility Xconn is compared in Fig. ^ with results from 
N t = 4 and our results from Nt = 6. This quantity is far less noisy than the disconnected 
contribution. The lack of a discernible peak in the Nt = 12 results is presumably another 
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consequence of the greater distance from the critical point. 

E. Crossover location 

To summarize, the observables we have considered place the crossover at N t = 12 in the 
range Q/g 2 = 5.65 — 5.70 at am q = 0.008 and in the range 5.75 — 5.80 at am q = 0.016, which 
we translate to 143 — 154 MeV p and 142 — 150 MeV p , respectively. This determination is 
compared with results for a wide variety of two-flavor simulations with both Wilson and 
staggered fermions in Fig | |9|JT3|]. Here we plot the temperature in units of the p meson 
mass vs the square of the ratio of the ir to p mass. For the staggered fermion simulations, we 
use the lightest local (non-Goldstone) pion. This state is expected to be degenerate with the 
Goldstone pion in the continuum limit, so progress toward the origin in this ratio measures 
restoration of flavor as well as chiral symmetry. Except for a few N t = 4 points from Wilson 
fermions, the consistency of these results (within about ten percent) is quite striking. 

III. SCALING TESTS 

A. Constituent Quark Free Energy 

The wide range of N t values now available makes possible an amusing analysis of the 
Polyakov loop variable, which measures the change in the free energy of the thermal ensemble 
due to the introduction of a point spinless test quark. This free energy difference, 

/(T,m g ) = -Tlog(ReP/3), (3.1) 

a function of the temperature T and light quark mass m q , includes the lattice-regulated 
ultraviolet-divergent self-energy of the point source, proportional to the inverse lattice spac- 
ing 1/a, and the free energy of the screening cloud of light antiquarks, quarks, and gluons, 
which we might call the "constituent quark free energy". Computing the self energy to 
leading order in perturbation theory, we have 
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f(T, m q ) = 2nC F a v -f/ a + f cq (T, m q ) . (3.2) 

where Cp = 4/3 is the color Casimir factor for the triplet representation, ay is the color 
fine structure constant appropriate for heavy quark bound states at the same lattice scale, 
and 

7 = iVf ^ 6 - 2 Ej=i cos(27r £; M /iV s ) + (Ma) 2 (3 ' 3) 
is the dimensionless static lattice propagator for a Debye-screened electrostatic gluon field 
evaluated at zero separation. Screening provides an arbitrary infrared cutoff. Although 
the ultraviolet divergent contribution 1/a is uniquely determined in the limit a — > 0, an 
infrared cutoff, here embodied in the Debye mass M, determines where the contribution 
from the point quark ends and the contribution from the screening cloud begins. Thus this 
approach offers no unique definition of a "constituent quark free energy" . Instead, within the 
framework of any consistent definition, Eq. ( |3.2| ) permits a separation of two contributions, 
one varying in a known way with the lattice scale 1/a and the other, unknown, but scale- 
invariant (at least in the continuum limit). Herein lies its predictive power: no matter what 
choice is made for the infrared cutoff, a measurement of Re P at two different values of N t 
can be used to predict ReP at any other N t . 

For present purposes, using N t = 1/aT, we choose a simpler approximate form 

f cq (T, m q ) = -T (log(Re P/3) + cN t ) (3.4) 

and adjust the dimensionless constant c by eye to achieve the rough scaling agreement 
shown in Fig. For each data set only values for the lightest available quark mass are 
used. Although the quark mass values m q /T are not the same from one Nt to the next in 
this figure, they are small (m q /T < 0.1) and the small variation would be expected to make 
little difference (of the order 10 MeV) in the free energy. The best value for c appears to be 
about 0.4, compared with values ranging from 0.26 to 0.35 expected from the lowest order 
perturbative self energy with a screening mass M = 3.2T. It is curious that at the crossover, 
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the free energy drops by about the 300 MeV expected in a constituent quark model with 
deconfinement at high temperature. 



B. Critical Behavior 

An SU(N) xSU(N) sigma model has often been proposed as the paradigm for the high 
temperature phase transition in QCD ||. For two light flavors (N = 2) one appealing 
alternative, consistent with numerical results, is that chiral symmetry restoration proceeds 
through a second order phase transition at zero quark mass with attendant 0(4) critical 
behavior. Flavor symmetry breaking in the staggered fermion scheme preserves only an 
exact 0(2) subgroup, so one may expect 0(2) critical behavior on coarse lattices. Away from 
the critical point, where fluctuations are unimportant, spin systems often exhibit mean field 
behavior. Kocic and Kogut have even proposed that mean field behavior could extend up to 
the critical point, as a consequence of the composite nature of the scalar sigma model fields 
as they appear in QCD jn$|. While it remains to be established whether their arguments 



apply to QCD, we consider the mean field alternative a plausible first approximation, since 
we do not know the extent of the Ginzburg scaling region. Thus even if 0(4) critical behavior 
is ultimately found in QCD at sufficiently large volume, it is possible that results imitate 
mean field behavior for lattices of the size considered in this study and at temperatures and 
quark masses that are not sufficiently close to the critical point. 

Karsch [^UJ and Karsch and Laermann || have analyzed the crossover, starting from the 



Fisher scaling hypothesis for the scaling of the critical contribution to the free energy 

f CIit (t,h) =b- d f CTit (byn,b yh h). (3.5) 

From this scaling behavior one derives a scaling relation for the critical contribution to the 
magnetization s = —df CT \ t /dh: 

s{t,h) = h l ' s y{x) (3.6) 
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where x = th^ 1 ^ 5 and y(x) is a scaling function. In QCD the quark mass plays the role 
of the magnetic field and (ijrtp^, the magnetization. Specifically, Karsch suggests using h = 
m q /T = am q N t and t = 6/g 2 — 6/g 2 (0, N t ), where g c (0, N t ) is the critical gauge coupling at 
zero quark mass for a particular Nt J20|- This identification leads to a successful accounting 
for the crossover location ( "pseudocritical point") 6/g pc (m q ) 2 for iV t = 4 at small quark 
mass. 

We extend the Karsch and Laermann analysis to a wider range of lattice spacings and 
test the scaling relation ( |3.6| ) directly. To anticipate the need to explore critical behavior, 
while simultaneously approaching the continuum limit, we propose an identification that 
avoids quantities with anomalous dimensions, but entails a translation to physical units: 

h = ml(m g , T = 0)/m%m q , T = 0) (3.7) 
*=[T-T c (0)]/T c (0) (3.8) 
s = h- 1 m q (W(m q ,T))/T i (3.9) 

The scaling relation (|3.6|) then gives a universal function 

y {x) = h- l -^ s m q (^(m q , T))/T 4 (3.10) 

with x = th~ l l@ 5 . The extra factor h~ l is needed to compensate for the quark mass factor 
m q . Let us apply this analysis to data for (vpip^ from several groups ||21|| . We begin by 
applying the following cuts to the world data: 

• m q /T = am q N t < 0.2 

• h < 0.5 

• -1 < x < 1 

• -0.5 < t < 0.5 

The first cut keeps parameters within range of the spectrum interpolation. The last three 
are an attempt to keep values within the region dominated by the scaling part of the free 

13 



energy. Since we have no a priori knowledge of the extent of this region, these choices are 
arbitrary. Indeed the cut on h is probably too generous, as we shall see, but we have set it 
high initially so that all of the N t = 12 data are included. 

We first carry out a mean field analysis of the data. The mean field magnetic equation 
of state y(x) is determined from the global minimum of the usual quartic free energy 

F(y,x) = by 4 + axy 2 — y (3-H) 

where a and b are adjustable scale parameters. The third adjustable parameter is the critical 
temperature T c (0) at zero quark mass, which we vary in 10 MeV increments. Our exploratory 
strategy is to seek the best agreement between world data (after the cuts listed above) and 
the resulting magnetic equation of state y(x), giving higher priority to data with higher N t . 
Setting the critical exponents to their mean field values, 5 = 1/3 and (3 = 1/2, and adjusting 
T c (0) to get the best qualitative agreement according to these preferences, we get the result 



shown in Fig. [10] with the choice T c = 140. Since the N t = 4 data appear to be outliers, it 
was necessary to omit them from the fit to the mean field magnetic equation of state. Apart 
from the N t = 12 data at higher quark mass, the agreement is surprisingly good. 

Our higher mass N t = 12 data are strikingly inconsistent with scaling. Indeed, as we 
have remarked in our analysis of the disconnected susceptibility, the am q = 0.016 data are 
taken at values of h higher than those of the rest of the sample included in the graph. Thus 
it is plausible that these data are taken too far from the critical region. 

An alternative strategy, the one used in |J, is to emphasize agreement among world data 
for the full range of N t without regard to the mean field magnetic equation of state. Our 



best result is shown in Fig. [TT| with the choice T c = 150. The mean field curve, based on 
fitting only over the range —0.25 < x < 0.25 clearly does not match the data as well as it 
did with T c = 140. 

The corresponding exercise with 0(4) critical exponents 5 = 4.82 and (3 = 0.38 is shown 
in Figs. O and [TBI. In this case the magnetic equation of state is not known in closed form. 



Instead, it was obtained from a simulation of an 0(4) spin model p2 |, using critical exponents 
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and the critical coupling reported by Kanaya and Kaya ||23|| . The level of agreement is 
comparable to that of the mean field analysis. 

At the present level of precision it is not possible to distinguish 0(4) from mean field 
critical behavior, let alone from 0(2) (not shown). Obviously a host of systematic errors, 
including finite volume effects and deviations from continuum scaling enter the analysis, so 
refinements are certainly needed before the method can serve as a definitive test of critical 
behavior. 



IV. SUMMARY AND CONCLUSIONS 

We presented results of a numerical simulation of QCD with conventional staggered 
fermions at nonzero temperature at the smallest lattice spacing attempted to date — 
approximately 0.1 fm. Numerical results are collected in Tables | - [TV] . Most of our data 
is obtained in simulations at N t = 12, but some new results are also presented for N t = 6. 
Locating the high temperature crossover at 0.1 fm for h = m^/m 2 p > 0.2 is difficult in the 
conventional observables (ReP) and {jpip). However, the baryon susceptibility and induced 
quark number show an abrupt rise consistent with what is found at lower N t . We use these 
signals to locate the crossover. At the lighter quark mass it is found to be in the range 
143 — 154 MeVp, consistent with results at lower N t . Judged in terms of the ratio h, our 
results at higher quark mass am q = 0.016 are farther from the critical point than the world 
data at smaller Nt. Thus the expected crossover peak in the disconnected chiral suscepti- 
bility is barely discernible in our am q = 0.016 data set, but it is present in our am q = 0.008 
data. We demonstrated a rough scaling analysis of the Polyakov loop variable in terms of 
a constituent quark free energy. Finally, (^pip^j shows remarkable agreement among world 
data with appropriate cuts over the range Nt = 6, 8, and 12 with a mean field magnetic 
equation of state and an assumed crossover temperature of approximately 140 MeV at zero 
quark mass. Our analysis is not refined enough to distinguish among 0(2), 0(4), and mean 
field scaling, however. 
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APPENDIX A: CONVERSION FROM LATTICE TO PHYSICAL UNITS 

The lattice scale is converted to physical units using a = am p / 770 MeV -1 following a 
method similar to that of Ref. plj. The numerator am p is found through a combination 
of interpolation and extrapolation from zero temperature spectral measurements described 
below. Of course, taking the physical value to be 770 MeV ignores shifts in this value as the 
quark mass is varied. We also require masses of the Goldstone pion and the second local 
pion m vr 2, which becomes degenerate with the pion in the continuum limit. Unfortunately, 
the highest value of 6/g 2 for which spectral data are available is 5.7, but our Nt = 12 data 
span the range 5.65 to 5.85, requiring extrapolation. For data at lower N t an interpolation 



suffices. Figure [14] indicates the parameter range of spectral data used for interpolation 



25fl . In selecting the world data sample, we have not attempted to correct for finite volume 
effects. 

For the interpolation we found a polynomial/spline expansion convenient: 

m^/rrip = R n/p (6/g 2 , am q ) = (am q ) 1/2 S n (6/ ' g 2 ) + am q T n (6/g 2 ) + (aniqf^U^Q/g 2 ) (Al) 
am p = R p (6/g 2 , am q ) = S p {Q/g 2 ) + am g T p (6/g 2 ) + (am q ) 2 U p (6/ ' g 2 ) (A2) 
am w2 = R^/g 2 , am q ) = S n2 (6/g 2 ) + amqT^jg 2 ) + (am q ) 2 U w2 {G/g 2 ) (A3) 
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(A4) 

where the S, T, and C/'s are natural cubic splines with three knots, 6/g 2 = 5.3,5.5,5.7 
selected to span the range of available spectral data. 

For the N t = 12 data over the range 6/g 2 = [5.65, 5.85], we chose an extrapolation based 
on two-loop asymptotic scaling: 

aA(6/g 2 ) = [8vr 2 / (3b g 2 )] b ^ b o exp([-87r 2 /(b g 2 )] (A5) 

where bo = 11 — 2Nf/3 and b\ — 51 — 19iV//3 and Nf = 2. In place of the bare lattice 
coupling, we used the tadpole improved value proposed by Mackenzie and Lepage || and 
generalized to two flavors by Bitar et al. |2f| : 

Q/gts = MO) = -2/log(D/3) - 1.99/tt (A6) 

where □ is the plaquette variable (normalized to 3 for the vacuum). Applying the extrap- 
olation requires knowledge of the plaquette over the range of extrapolation. Since zero 
temperature values are unavailable over this range, we used the plaquette values measured 
in our high temperature study, using a spline interpolation formula in the same pattern as 
the masses above, but with knots at 6/g 2 = 5.65,5.75,5.85. At our small lattice spacing, 
the temperature effect on the plaquette is negligible for this purpose. To avoid introducing 
distortions in our predicted mass values resulting from a switch in method from the spline 
interpolation to the extrapolation procedure, we elected to use only an extrapolation, based 
on interpolated spectral data at 6/g 2 = 5.65, the lowest value in our data set, rather than 
at 6/g 2 = 5.70, the highest parameter for which spectral data is available. One must allow 
for a rescaling of the quark mass in the extrapolation. Thus to determine the starting point 
for extrapolation, we first solve for the scale parameter s 

s = A{f3 cS [a(6/g 2 ,am q )]} 

A{/3eff[n(5.65,am ff /s)]}' 1 ' 

The extrapolated masses are then obtained by rescaling from the interpolation at 6/g 2 = 
5.65: 
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m n /m p = R n / p (5. 65, am q /s) (A8) 
arrip = sR p (5.Q5,am q /s) (A9) 
am^ = sR n 2 (5.65, am q /s). (A10) 

The measured mass values at 5.70 can then be used as a check of the extrapolation. In 
Fig. |15| we compare the empirical formula for hadron masses with the measured values. The 
extrapolation from 6/g 2 = 5.65 predicts slightly higher mass values at Q/g 2 = 5.7 than 
measured. The agreement is certainly adequate for present purposes. To refine the method 



would require a more detailed study of finite volume effects |27 . 
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FIGURES 
FIG. 1. Chiral condensate (dip) vs 6/g 2 . 

FIG. 2. Polyakov loop vs 6/g 2 . 

FIG. 3. Fuzzy Polyakov loop vs 6/g 2 . 

FIG. 4. Nonsinglet baryon susceptibility vs 6/g 2 for Nt = 4,6,8,12. Free lattice quark values 
for each N are indicated by horizontal lines on the right of the plot. 

FIG. 5. Induced quark number. Octagons are for fixed lattice spacing. 

FIG. 6. Disconnected chiral susceptibility. Crosses for Nt = 4, am q = 0.0375 ||, octagons for 
N t = 4, am q = 0.02 ||, squares for Nt = 6, am q = 0.0125, diamonds for N t = 12, am q = 0.016, 
and bursts for Nt = 12, am q = 0.008. 

FIG. 7. Connected chiral susceptibility. Same notation as Fig. ^. 
FIG. 8. Crossover temperature in units of the rho meson mass vs the square of the pi to rho 



mass ratio. For staggered fermions the mass of the lightest non-Goldstone pion is used [13|. The 
vertical line locates the physical ratio. 

FIG. 9. Constituent quark free energy as defined in text. 

FIG. 10. Scaled {ipvp} vs scaled temperature with mean field critical exponents at T c (0) = 140 
MeV. 

FIG. 11. Same but with T c (0) = 150 MeV. 
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FIG. 12. Same but with 0(4) critical exponents and T c (0) = 130 MeV. 

FIG. 13. Same but with 0(4) critical exponents and T c (0) = 140 MeV. 

FIG. 14. Range of parameter pairs (6/g 2 ,am q ) for two-fiavor staggered fermion thermodynamic 
simulations (various plot symbols as indicated) and spectral data (symbol "S") used to set the 
lattice scale. 

FIG. 15. Empirical fit to the ir, p, and 7T2 masses in the compilation of available spectral data. 
Open symbols (square, octagon, diamond) and error bars show the respective measured values. 
The cross, burst, and plus show the respective empirical estimates. The -K2 values at 5.5, 5.6, and 
5.7 are offset for clarity. 
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TABLES 

TABLE I. Summary of results for N t = 12, am, = 0.008: chiral condensate, Polyakov loop, 
fuzzy Polyakov loop and induced baryon charge. 



Q/g 2 (tpip) ReP ReF Q 

5.65 0.0393(3) 0.0033(5) 0.00051(6) -0.6(2) 

5.70 0.0311(5) 0.0094(15) 0.00155(11) -0.05(7) 

5.725 0.0290(2) 0.0118(17) 0.00178(15) -0.16(7) 

5.75 0.0267(2) 0.0168(13) 0.00255(12) -0.03(3) 

5.80 0.02424(10) 0.0212(10) 0.00359(10) -0.03(2) 

5.85 0.02262(7) 0.0254(12) 



TABLE II. Summary of results for Nt = 12, am q = 0.008: nonsinglet and singlet baryon 
susceptibility and disconnected and connected chiral susceptibility. 

6/s 2 Xns/T 2 Xs/T 2 Xdisc/T 2 Xco„n/T 2 

5.65 0.6(2) 0.4(5) 54(10) 480(8) 

5.70 1.40(13) 1.3(5) 80(20) 466(12) 

5.725 1.44(11) 1.5(4) 29(6) 469(9) 

5.75 1.52(10) 1.6(3) 25(7) 458(7) 

5.80 1.87(8) 1.8(4) 6.2(1.1) 435(6) 
5.85 
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TABLE III. Summary of results for N t = 12, am q = 0.016: chiral condensate, Polyakov loop, 
fuzzy Polyakov loop and induced baryon charge. 

6/g 2 ReP ReF Q 

5.65 0.0672(2) 0.0034(4) 0.00050(5) -0.2(2) 

5.70 0.0599(2) 0.0046(9) 0.00065(8) -0.4(2) 

5.75 0.0522(3) 0.0136(8) 0.00205(9) -0.12(4) 

5.80 0.0491(3) 0.0164(12) 0.00225(15) -0.10(3) 

5.85 0.04516(14) 0.0236(13) 0.00368(14) 0.01(3) 



TABLE IV. Summary of results for N t = 12, am q = 0.016: nonsinglet and singlet baryon 
susceptibility and disconnected and connected chiral susceptibility. 

6/ 5 2 Xns/T 2 Xs/T 2 Xdisc/T 2 Xco„n/T 2 

5.65 0.48(15) 0.9(4) 31(6) 434(6) 

5.70 0.57(12) 0.4(5) 24(6) 427(5) 

5.75 1.14(12) 0.4(4) 26(3) 427(7) 

5.80 1.21(10) 1.3(3) 30(5) 404(4) 

5.85 1.58(11) 1.2(4) 8(2) 382(5) 
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